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0.  Introduction 

There  is  now  considerable  interest  in  initial -boundary  value 
problems  for  various  systems  of  partial  differential  equations  arising 
in  fluid  dynamics.  This  interest  stems,  primarily,  from  efforts  to 
create  useful  computational  models  of  various  processes  for  the  purposes 
of  prediction  (atmospheric  processes,  ocean  circulation,  etc.)  and  the 
detailed  study  of  various  phenomena  (convection,  flow  in  wind  tunnels, 
lee  waves,  eddies,  etc.).  Such  calculations  are  not  new.  As  these 
computational  models  have  become  more  accurate  difficulties  with  the 
boundary  conditions  have  become  more  evident.  This  has  led  first  to  the 
examination  of  the  various  discretizations  used  and  then  back  to  the 
differential  equations  whose  approximate  solutions  are  sought. 

Such  a backward  sequence  of  events  may  seem  surprising.  Naturally, 
the  initial -boundary  value  problems  for  the  differential  equations  should 
have  been  carefully  examined  first  since  we  cannot  expect  our  approximations 
to  be  reasonable  if  they  approximate  a problem  which  does  not  have 
reasonable  solutions.  The  reason  it  has  gone  this  way  is  clear. 
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It  is  natural  to  first  examine  the  evidence  where  it  appears  and,  as  usual, 
the  computations  have  been  ahead  of  the  analysis.  The  initial-boundary 
value  problems  for  these  systems  of  differential  equations  are  not  easy 
to  analyze;  and,  in  fact,  adequate  tools  for  a rather  complete  analysis 
have  only  recently  become  available  stemming  from  the  work  of  Kreiss  [12,13] 

The  current  interest  has  resulted  in  several  works  based  on  the 
classical  energy  method  (e.g.,  Elvius  and  Sundstrdm  [9],  Davies  [5 ,6], 
de  Rivas  [ 7 1 and  Dutton  [8J),  which  follow  the  earlier  work  of  Serrin 
[17],  Sundstrbm  [19]  and  Campbell  [1].  However,  this  method  only  works 
for  a limited  class  of  equations  and  boundary  conditions.  Some  authors 
have,  unfortunately,  made  unallowable  assumptions  (over-specification  of 
boundary  conditions,  omission  of  terms,  etc.)  in  futile  attempts  to  make 
their  problems  fit  into  this  class.  We  will  discuss  some  instances  of 
this  in  detail.  This  seems  to  be  a real  hazard  in  the  use  of  the  energy 
method  since  the  effects  of  such  assumptions  are  often  well-buried  in 
intermediate  estimates  and  consequently  overlooked. 

We  begin  with  a general  discussion  of  well-posedness  followed  by 
a review  of  properties  of  the  adiabatic,  inviscid  Eulerian  equations  of 
fluid  dynamics  (system  A).  We  then  study  two  approximations  of  the 
Eulerian  equations:  the  hydrostatic  "primitive  equations"  of  meteorology 

(system  Bl)  and  the  incompressible,  anelastic  equations  (system  32). 

Finally,  we  discuss  the  shallow-water,  or  barotropic,  equations  (system  C) 
which  can  be  considered  as  a further  simplification  of  system  Bl  or  B2. 

It  is  interesting  to  consider  these  equations  in  this  order  so  that  the 
effect  of  each  successive  approximation  can  be  observed.  The  systems  A 
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and  C are  symmetrizable , hyperbolic  systems  but  systems  B1  and  B2  are 
not  hyperbolic.  These  facts  have  profound  influence  on  the  well-posedness 
of  initial  boundary-value  problems  for  these  systems. 

We  consider  two  types  of  boundary  conditions  which  arise 

naturally  in  many  situations.  Most  of  our  analysis  will  deal  with 

certain  quarter-space  problems  but  we  will  always  have  the  following 

2 

underlying  situation.  Let  fi  c R be  an  open,  connected  region  with 
smooth  boundary,  dfi,  and  (i  = 0 U dn.  We  will  consider  the  system  C 
on  the  domain  = Q x [0,T]  and  the  systems  A,  B1  and  B2  on  the  domain 
= 0 * I * [0,T]  where  I = [0,<»)  or  I = [0,1].  The  two  types 
of  boundary  conditions  we  consider  on  dfi  are:  (l)  rigid  wall  boundaries 

and  (2)  open  boundaries.  The  rigid  wall  case  corresponds  to  a physical 
situation  which  requires  the  normal  velocity  to  vanish  at  the  boundary 
and  is  the  simpler  of  the  two  types.  This  situation  is  often  encountered 
in  oceanography.  Open  boundaries  occur  in  limited  area  forecasting, 
wind  tunnel  flow,  and  studies  of  small  scale  or  local  phenomena  in 
meteorology  and  oceanography.  Open  boundaries  do  not  arise  from  a 
natural  physical  situation  and  a suitable  form  for  the  boundary  con- 
ditions is  not  obvious.  Boundary  conditions  which  do  not  introduce 
boundary  layer  phenomena  are  usually  wanted  in  this  case.  That  is,  these 
boundary  conditions  should  determine  the  interior  flow  as  though,  in  fact, 
the  boundaries  were  not  there  at  all.  In  each  case  we  give  necessary 
conditions  for  the  form  of  the  boundary  conditions  in  order  that  the 

problems  be  well-posed.  We  will  also  give  particular  boundary  conditions 
which  yield  well-posed  problems.  We  show  that  systems  A,  B2,  and  C can 
be  treated  satisfactorily  and  that  system  B1  is  ill-posed  for  local, 

pointwise  boundary  conditions.  For  the  linearized, constant  coefficient 
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versions  of  systems  B1  and  B2,  well-posed  boundary  conditions  are 
given.  It  seems  reasonable  to  conjecture  that  these  boundary 
conditions  also  yield  a well-posed  problem  for  system  B2.  However, 
the  corresponding  boundary  conditions  for  system  B1  cannot  be  easily 
implemented  for  the  general  problem. 

As  already  mentioned,  many  of  the  presently  used  boundary  con- 
ditions specify  more  data  than  is  allowed.  These  specifications  preclude 
the  existence  of  smooth  solutions  except  in  very  special,  unrealistic 
situations  where  the  exact  solution  is  known  on  the  boundary  without  error. 
Errors  must,  however,  be  expected  in  the  boundary  data  arising  from  errors 
in  measurement,  the  use  of  constant  boundary  values,  or  from  computations 
over  larger  regions  if  some  telescoping  grid  technique  is  used.  We 
will  discuss  the  implications  of  such  overspecifications. 

Most  of  our  analysis  deals  with  inviscid  systems  of  equations. 
Viscous  terms  are  added  to  these  equations  in  many  forecast  integrations. 
They  are  often  motivated  physically  as  representing  "eddy  diffusion" 
of  momentum  and  potential  temperature.  The  effect  of  these  terms  on  the 
main  part  of  the  solution  is  usually  small.  The  real  motive  for  including 
them  is  often  non-physical.  Since  the  equations  are  nonlinear, 


initial  long-wave  phenomena  can  produce  shorter  wave  phenomena  which  cannot 
be  accurately  represented  by  the  approximation  used.  To  prevent  result- 
ing aliasing  errors  and  nonlinear  instability  the  computational  method 
should  be  provided  with  a dissipative  filter  term,  and  the  most  primitive 
form  of  filter  is  just  such  an  "eddy  diffusion"  term. 

In  both  cases,  the  viscous  coefficients  are  so  small  that  we 
should  expect  the  boundary  conditions  to  be  close  to  those  valid  for 
the  corresponding  inviscid  system.  The  viscous  equations  do,  however, 
require  additional  boundary  conditions,  and,  as  an  effect,  viscous 

k _ 


boundary  layers  may  occur  at  the  boundaries.  Such  boundary  layers  may 
sometimes  be  appropriate,  as  in  the  rigid  wall  situation.  However, 
at  open  boundaries,  they  are  inappropriate.  We  shall  therefore  discuss 
the  formulation  of  boundary  conditions  when  viscous  terms  are  included 
and  show  how  these  conditions  can  be  chosen  so  that  no  singular  boundary 
layers  result  as  the  cofficients  of  the  viscous  terms  tend  to  zero. 


1.  Well-Posedness 

Our  main  goal  is  to  establish  the  existence  or  non-existence  of 
certain  a priori  estimates,  or  energy  inequalities,  valid  for  the  solutions 
of  the  various  initial- boundary  value  problems  under  consideration  here. 

In  this  section  we  discuss  the  form  of  these  estimates  and  some  implica- 
tions that  follow  from  them. 

For  the  purposes  of  this  discussion  let  us  write  our  problems 
in  the  general  form 


ll 

& 

in 

fi  x 

I x 

[0,T] 

(the  differential  equation) 

, A = Ho 

in 

n x 

I 

(the  initial  conditions ) 

t=o 

u1  = Su11  + g 
^ 

in 

dfi  x 

I x 

[0,T] 

(the  boundary  conditions ) 

where  L is  a partial  differential  operator;  u,  F and  u^  are  vector 
functions  of  dimension  k,  u = (u^u11)^;  u1  and  g,  are  of  dimension  i, 
u is  of  dimensionsion  k -l'}  and  S is  a real  l * (k-£)  matrix. 

For  a linear,  first  order  hyperbolic  equation  in  R , -L  would  take 
the  form 


(1.2) 


Aj (*»*) 


B(x,t)  . 
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The  boundary  conditions  express  the  components  u1  of  u in  terms 
of  the  components  u11  of  u and  the  given  function  g.  The  matrix  S 
can  be  thought  of  as  a generalized  reflection  operator.  Changes  of 
variable  may  be  necessary  to  bring  certain  desired  boundary  conditions 
into  the  form  given  here  but  this  does  not  introduce  any  essential 
restriction.  The  partitioning  u = (u1,  u^)t  is  based  upon  the 
characteristic  variables,  or  Riemann  invariants,  of  the  problem,  that  is,  the 
components  of  u^  can  be  called  "incoming"  quantities  and  those  of 
u^  can  be  called  "outgoing"  quantities.  This  partitioning  of  u 
will  be  discussed  in  detail  for  the  problems  we  consider. 

The  estimates  we  seek  are  of  the  form 

IX  [0,T]  + 5^SfiXIX[o,T]  + 

KT 

1 e dl£  Hnx  ix  [o,t]  + %^dnx  ix  [ o,t]  + ^Hfix  i^ 

2 2 

where  the  norms  are  L norms  or  weighted  L norms  over  the  regions  indicated 
by  their  subscripts,  K 0 is  a constant  independent  of  T,  and  & = 0 
or  1.  We  will  refer  to  (1.3)  with  3 = 0 as  the  weak  form  of  (1.3).  The 
differences  in  the  properties  of  solutions  u which  satisfy  only  the 
weak  form  (1.3)  from  those  which  satisfy  (1.3)  with  & - 1 are  discussed 
by  Kreiss  [133*  We  will  not  belabor  the  distinction  here  and  be  satisfied 
with  the 
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Definition  1.1.  We  will  say  that  the  problem  (l.l)  is  well-posed  if 
the  estimate  (1.3)  holds  for  all  solutions  u of  (l.l)  with  F,  ir. 
and  g in  L2. 

Uniqueness  and  stability  with  respect  to  perturbations  in  the 

data  follow  from  the  estimate  (l>3).  We  refer  to  the  works  of  Kreiss  [12, 

13],  Majda  and  0sher[l4]and  Strikwerda  [ 18]  for  detailed  discussions 

0 

of  the  particular  weighted  L norms  and  the  general  theory  for  systems 

of  hyperbolic  and  incompletely  parabolic  equations. 

The  equations  we  are  considering  are  all  quasi-linear. 

However,  we  can  obtain  our  a priori  estimates  a posteriori  over  intervals 

[0,T]  when  a smooth  solution  exists,  i.e.,  we  can  consider  coefficients 

A (u(x,t),x,t)  as  functions  of  x,t  if  u(x,t)  is  known.  Furthermore, 

0 

iterations  based  upon  the  linearized  variational  form  of  the  problems 
can  be  used  to  establish  existence  for  those  t-intervals  where  the 
iterations  converge.  We  will  not  pursue  this  here,  but  rather  assume 
the  existence  of  smooth  solutions  over  the  interval  of  consideration. 

We  must  prescribe  boundary  conditions  that  do  not  preclude  the 
existence  of  such  smooth  solutions.  This  is  the  case  if  too  many 
conditions  are  specified.  Too  few  conditions  preclude  uniqueness,  of 
course. 

The  results  by  Kreiss  [12]  and  Strikwerda  [18]  also  show  that 
the  systems  (a)  and  (c)  are  stable  to  perturbations  by  lower  order  terms. 
Th^s  implies  that  we  need  not  consider  the  effects  of  terms  such  as 
undifferentiated  frictional  terms  and  coriolis  forces  in  our  analysis. 
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It  is  essentially  due  to  this  fact  that  the  analysis  of  variable  coefficient 
problems  can  be  reduced  to  that  of  corresponding  constant  coefficient 
problems  via  the  construction  of  appropriate  pseudo-differential  operators 
iKreiss  f 12  ],  Taylor  [20],  Majda  and  Osher  [14] , Strikwerda  [lB]).  This 
stability  property  also  allows  us  to  reduce  problems  on  our  general 
domain  £2  with  smooth  boundaries  to  families  of  quarter-plane  problems 
by  making  local  changes  of  coordinates  such  that,  e.g.,  bft  is  mapped 
into  x^  = 0 and  fi  into  > 0.  Such  mappings  only  introduce  new 
terms  which  are  of  lower  order.  More  detail  about  this  process  can  be 
found  in  Majda  and  Osher  [l4]  and  Strikwerda  [13]. 

Existing  theoretical  results  cover  problems  with  a smooth  non- 
characteristic boundary  for  classes  of  equations  which  include  A and 
C and  their  modifications  resulting  from  the  inclusion  of  the  usual  eddy 
viscosity  terms  (Kreiss  [12],  Majda  and  Osher  [l4],  Strikwerda  [18]). 
Extensions  to  problems  in  regions  with  corners  and  uniformly  characteristic 
boundaries  have  been  studied  by  Majda  and  Osher  [14].  However,  the 
important  case  where  the  velocities  change  sign  on  the  boundary  and  do 
not  vanish  in  a neighborhood  of  such  a boundary  point  Is  not  covered  by 
existing  theory.  This  often  occurs  in  the  applications  we  consider, 
e.g.,  the  solid-wall  type  of  boundary  conditions, and  when  the  flow  direc- 
tion reverses  to  change  an  inflow  or  outflow  section  of  the  boundary 
to  an  outflow  or  inflow  section,  respectively.  There  must  be  character- 
istic points  on  any  smooth  boundary  of  a simple  connected  region  with 
open  boundaries  which  has  both  inflow  and  outflow  sections  of  the  boundary. 
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We  cannot  treat  the  influence  of  such  points  on  dfi  here  but  conjecture 
that  no  important  modifications  are  usually  necessary  for  problems  like 
those  we  treat  here. 

We  will  use  both  the  classical  energy  method  and  Kreiss'  normal 
mode  analysis  to  establish  the  well-posedness  of  these  problems.  The 
solid-wall  boundary  problems  are  all  treated  using  the  energy  method  ^ 
which  provides  us  with  estimates  of  the  form  (1.3)  directly.  Some 
boundary  conditions  for  the  open  boundary  problems  can  be  treated  in 

this  way,  but,  in  general,  we  must  use  normal  mode  analysis  for  these 
problems . 


2.  The  Eulerian  Questions  (System  A). 

The  basic  hydrodynamic  and  thermodyanmic  laws  governing  the 
motion  of  an  adiabatic  and  inviscid  fluid  are  given  by  the  Eulerian 
equations 


(2.l) 


^u  + af?p  + F = 0 

^ a - of?* u = 0 
It  ? + = 0 


where  u is  the  three-dimensional  velocity  vector,  u = (u1,Ug,u^)t, 

Of  is  the  specific  volume,  and  p the  pressure  of  the  fluid;  y = Cp/Cv 
is  the  lapse  rate  of  the  fluid,  F represents  zero-order  and  forcing 
terms,  e.g.,  coriolis  and  gravity  forces,  and 
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d _ d „ d y d 

dt  _ dt  + U dt  , , Uj  dx 

J **•  J 


In  vector  notation,  the  equations  are 


(2.3) 


KT  q + t A (q)  q + Fig)  = 0 


* * ' A “j  ^ ^ * 


where  q = (u^u^u  ,at,p)  , and 
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The 


matrices  A.  all  have  real  eigenvalues  u , u.,  u.,  u.  + c,  and 
" J J j J 

u^-c,  with  distinct  eigenvectors,  c = (pra)1'2  is  the  sound 


speed  of  the  fluid.  The  matrices  are  not  symmetric  but  it  is  easy  to 
find  a symmetric,  positive  definite  matrix 


1 
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R = a 


10  0 
0 i 0 


0 0 1 


oooa  pr/a 


0 0 0 


(1+a2)  -SL 

pr 


*-l  -1 

= T T 


where  a is  an  arbitrary  real,  non-zero  parameter  such  that  the  trans- 
formed matrices  T-1A.T  are  all  symmetric.  The  system  A is  thus  a 

J 

quasi-linear  system  of  hyperbolic  partial  differential  equations,  see 
e.g. , Courant  and  Hilbert  [ 4 ].  Since  no  closed-form  expression  for 
the  solution  to  this  system  is  known,  a rigorous  evaluation  of  the 
effects  of  different  approximations , inhomogeneous  terms,  and  boundary 
conditions  is,  in  general,  impossible.  For  the  problems  we  are  con- 
sidering the  solutions  are  usually  continuous  and  smooth.  The  matrices 
A j and  R are  then  also  smooth.  As  long  as  the  deviations  g>'  (x,t) 
from  the  exact  solution  £(x,t)  are  small,  they  should  then  approximately 
satisfy  the  linearized  variational  equations 


(2.4) 


Its'  + Ji1V*)  &l'  + cl'+l'  = 0 


where  a*  = (uj,  u£,  ul,  a’,  p* )t, 
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and  the  term  F'  can  represent  other  low-order  terms.  This  system  is 
a linear  hyperbolic  system  in  q' , and  the  well-posedness  can  be  studied 
by  either  the  classical  energy  method  or  by  Kreiss'  normal  mode  analysis. 
In  the  energy  method,  the  basic  idea  is  to  show  that  a suitable  norm  for 
satisfies  a growth  equation  of  the  form 


(2.5) 


£ Ik*  II  < Klk'll  + IM 

vjw  />*  f>*  r** 


where  ||*||  is  an  inner-product  norm  |k'  ||  = (fQ  q '*  Mq'  dx1dx2dx^ ) 
equivalent  to  the  Euclidean  L2-norm  |q'|2  dx..dx^dx,  )^2.  We  can 

u < ± 2 j 

show  the  well-posedness  of  the  pure  initial-value  or  Cauchy  problem 
2 

in  the  L -norm  for  the  Eulerian  equations  by  choosing  M as  the  matrix 
R given  above.  From 
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(2.6)  (q'  Rq'  ) 


£ U'Xa' 

3=1  ° j J~ 


*RA  .q' ) + q'*  + £ RA.-RC-C*R}  q' -q' *RF' -F' *Rq' 


A j' 


we  get 


(2*7) 

= f[%'*  (If  + ^ 57"  GO-RC-CR}  q'  - q'*RF‘-F*Rq,]ta1dx2dx 

J ,i=l  j J ~ ~ ~ ~ ~ 


and  since  R,  RA.,  etc.,  are  bounded,  slowly  varying  matrices,  we  can 
J 

easily  establish  a growth  equation  of  type  (2.5  )• 

This  inequality  and  energy  norm  is  essentially  equivalent  to  that 
used  by  Serrin  [17]  in  his  uniqueness  proof  for  compressible  fluids. 

He  wanted  the  estimates  to  be  valid  for  large  deviations  , which 
complicated  the  structure  of  the  proof.  However,  at  one  step  in  his 
calculations,  (eq.  25),  he  had  to  make  an  assumption  which  essentially 
meant  that  the  deviations  must  be  small.  Furthermore,  in  the  analysis  of  the 
limited-area  case,  an  over-specified  set  of  boundary  conditions  was 
used  on  the  inflow  portion  of  the  boundary,  thereby  invalidating  this 
part  of  the  proof.  The  computational  effects  of  such  an  overspecification 
will  be  discussed  later. 

For  the  initial-boundary  value  problem,  we  can  use  the  growth 
equation  (2.6)  as  before  if  and  only  if  the  boundary  integral  provides  a 
non-negative  contribution,  i.e.,  if 
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(2.8) 


RA  q’  ds  > 0 

n A — 


where  A 


= Z^  ^ A nt  • e . , n i 
n .1=1  ~ 


is  the  unit  vector  in  the  outward  normal 


direction,  and  e.  is  the  unit  vector  in  the  x .-direction.  Note  that  if 

J 

the  boundary  conditions  are  such  that  this  inequality  is  satisfied,  then 
well-posedLness  is  proved,  otherwise  no  conclusions  can  be  made. 

The  integrand  q'  RA  q'  is  a quadratic  form  in  the  five  variables 
UiJ  U2 ’ a’  > 8111(1  P' • However,  the  number  of  boundary  conditions  is 
only  equal  to  the  number  of  inward  characteristics,  that  is,  the  number 
of  negative  eigenvalues  of  A , see  Kreiss  [12].  These  boundary  conditions 
must  be  such  that  the  related  combinations  of  £'  (the  Riemann  invariants) 
are  given  in  terms  of  known  quantities  and  combinations  corresponding  to 
outward  characteristics. 

The  initial-boundary  value  problems  for  the  Eulerian  equations 
arise  from  two  different  situations  which  must  be  studied  separately: 

1.  A solid-wall  boundary.  Here,  the  physical  boundary  condition  is 

that  the  normal  velocity  v = YP.  , u . n^  • e . should  vanish  at  the 

n j=l  ~j 

boundary.  This  condition  is  consistent  with  the  number  of  inward 
characteristics  (one).  Since  it  also  gives  q'  RA  q'  = 0 identically, 
the  well-posedness  of  the  initial-boundary  value  problem  follows 
directly. 

2.  An  "open  boundary,"  or  a boundary  located  in  the  interior  of  a body 

of  fluid.  In  this  case  the  normal  velocity  is  non-zero  on  the 

boundary,  except  at  certain  points.  For  supersonic  fluids,  all 

quantities  u^,  u u^,  a,  and  p should  be  specified  at  inflow 
* 

points,  giving  q'  RA  q'  = 0.  At  outflow  points,  no  boundary 


r 


conditions  should  be  prescribed  and  q'  RA  q'  is  always  > 0. 

For  subsonic  fluids , has  four  negative  eigenvalues  with  distinct 

eigenvectors  at  inflow  points  so  that  four  quantities  have  to  be 

specified.  At  outflow  points,  the  number  of  negative  eigenvalues 

and  quantities  to  be  prescribed  is  only  one.  In  both  cases  no 

obvious  physical  boundary  conditions  are  known.  The  characteristic 

combinations  of  q'  may  be  found  from  the  eigenvectors  of  A , but 

n 

it  is  easier  to  proceed  as  follows. 


The  quadratic  form 


q'  RA  q' 


una  ^u^2  + uj_2  + uj|  + a2a,2pr/a+  (l+a2  )p'2a/pr  + 2a2a'p' } + 2uJ, 


(u|  1 and  u|p  are  the  velocity  components  in  two  orthogonal  tangential 
directions ) can  be  rearranged  as 

= una  1{u]^  + u'^+  a2[a' (pr/a)1^2  + p' (a/pr)1//2]2) 

+ \ + un ^ a 1(u,l  + P'  (a/pr )1^2)2  - \ (c-u  ) a”1(u'-p' (a/pr )1//2)2 


As  before  , c = (pra)  ' is  the  basic  sound  speed. 

At  an  inflow  part  of  the  boundary  (un  < 0),  the  characteristic 
combinations  corresponding  to  negative  eigenvalues  of  A^  are  thus 
u|  i > u^  - p' (a/pr)1^2)  and  a'Cpr/q)1^2-*-  p'Ca/pr)1^2  = cB'/p,  where 

6 is  the  "potential  temperature"  9 = (pl”"1//r/R)Qp1/r . The  four  boundary 


conditions  should  give  relations  of  the  form 


1 


= a^u^  + p’ (a/pr  )1//2)  , 
ui2  = a2^un  + p' (a/pr)1^2)  , 
- p ' (a/pr)1^2  = a3^un  + P '(a/pr)1/2), 
a' (pr/a)1^2  + p'(a/pr)1//2  = a^(u^  + p' (a/pr)1^2) 


for  deviations  from  the  basic  solution. 

Before  studying  specific  examples  of  boundary  conditions  giving 
relations  of  this  type,  it  should  be  noted  that  not  all  such  conditions 
give  well-posed  problems.  With  the  classical  energy  method,  we  can 
actually  only  prove  well-posedness  when 


(2  2,  2 2 \ . 1 / ,2  l , , 

u (a  + a0  + a a.  ) + — (c  - u ;a  ^ — (c  + u ) 

n 1 2 4 2 n 3 ~ 2 n 


since  only  then  is  q’  RA  q'  > 0 at  these  inflow  points.  In  this  expres- 
2 

sion  a is  the  arbitrary,  real  and  positive  parameter  of  the  matrix  R. 

If  we  want  to  investigate  the  well-posedness  for  other  values  of 
a^,  a2,  a^,  a^,  we  have  to  use  the  normal  mode  analysis  of  Kreiss, 
see  Oliger  and  Sundstrom  [15]. 

The  classical  energy  method  certainly  works  if  the  boundary 
conditions  are  such  that  a^^  = = a^  = a^  = 0.  We  cam  obtain  these 

relations  by  prescribing  (at  inflow  points)  the  two  tangential  velocity 
components,  the  potential  temperature  0,  and  the  combination 

un  • A (Prot)l/2‘ 

The  inequality  is  not  satisfied  if  we  have  a1  = a^  = a^  = 0 


and  a^  = + 1,  that  is,  if  we  try  to  give  the  tangential  velocity 


16 


components,  the  potential  temperature,  and  either  the  normal  velocity 
or  the  pressure  at  these  inflow  points.  A third  possibility  for  which 
the  energy  method  does  not  work  is  a^  = a^  = 0,  = -1,  and  a^  = 1, 

that  is,  if  we  give  all  three  velocity  components  and  the  specific 
volume  a.  As  shown  in  Oliger  and  Sundstram  [ 15]  using  the  normal 
mode  analysis  technique,  this  last  combination  actually  gives  a well-posed 
problem. 

At  the  outflow  part  of  the  boundary  (un  > 0),  only  one  quantity 
should  be  prescribed.  It  should  give  a relation  of  the  form 

u^  - p'  (a/pr  )1//2 

= 'b1u]1  + + b3^un  + p’  (a/pr)1/2)  + b^(a’  (pr/a)1/2  + p’  (a/pr)1/2) 

for  the  linearized  variational  equations.  It  is  simplest 
to  prescribe  the  normal  velocity  so  that  b1  = b£  = b^  = 0,  b^  = -1. 

The  well-posedness  of  this  boundary  condition  follows  immediately  from 

V- 

the  positivity  of  q*  RA  q* . We  may  instead  prescribe  the  pressure  p, 

~ n 

corresponding  to  b^  = bg  = b^  = 0,  b^  = 1,  or  actually  any  combination 

of  u and  p such  that  (u1  -p'  (a/pr  )2  (c-u  ) < (u’+p*  (a/pr^^J^c+u  ). 
n n n — n n 
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3-  Basic  Approximate  Forma  of  the  Eulerian  Equations  (Systems  B1  and  B2 ) 

Although  the  Eulerian  equations  are  the  fundamental  system  of 
equations  for  most  fluid  flow  problems,  they  are  often  used  in  modified 
approximate  forms.  In  fact,  they  have  almost  never  been  used  in  the 
complete  form  given  in  the  last  section  for  geophysical  calculations. 

The  reason  for  this  is  simply  a matter  of  economics  and  time,  which 
are  not  unrelated.  To  compute  an  accurate  approximate  solution  of  the 
Eulerian  equations  for  a relatively  small  problem  requires  quite  a lot 
of  computer  time.  These  equations  are  an  extremely  "stiff"  system  of 
hyperbolic  equations  with  a wide  range  of  eigenfrequencies  and  character- 
istic phase  velocities.  The  ratios  between  the  largest  and  smallest 

eigenvalues  of  the  coefficient  matrices  A.  in  (2-3 ) are  often  as 

J 

large  as  102  or  10  . The  high-frequency  eigensolutions  (sound-waves) 
are  often  absent  in  the  initial  data  and  the  solution, but  their  presence 
in  the  set  of  eigensolutions  imposes  a severe  upper  limit  for  the  time- 
step  in  explicit  numerical  integration  procedures.  Implicit  techniques 
that  do  not  suffer  from  this  difficulty  lead  to  data  structures  which 
are  difficult  to  manage  and  systems  of  nonlinear  equations  that  are 
expensive  to  solve. 

A second  special  aspect  of  many  geophysical  problems  is  the 
strong  balance  between  the  gravitational  and  vertical  pressure  gradient 
forces  which  is  responsible  for  the  basic  stratification  of  the 
atmosphere  and  oceans.  The  vertical  acceleration  terms  are  usually  much 


less  than  10  times  either  of  these  terms.  Even  if  we  first  subtract 
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the  time -independent  part  of  the  pressure  field,  we  must  still  know 
the  specific  volume  and  pressure  extremely  accurately  in  order  to  compute 
the  time  derivative  of  the  vertical  velocity  with  even  moderate  accuracy. 
Similar,  hut  less  extreme,  balances  exist  in  the  other  equations.  In 
the  two  remaining  equations  of  motion  there  is  a near  balance  between 
the  horizontal  pressure  gradient  and  coriolis  terms,  and  in  the  continuity 
equation  there  is  a near  balance  of  the  components  of  divergence.  These 
relations  are  often  summarized  in  the  statement  that  the  atmosphere  is, 
to  a large  extent,  not  only  quasi -hydrostatic  but  also  quasi -geostrophic 
and  quasi -nondivergent. 

To  obtain  a reasonable  computational  problem  we  must  either: 

(l)  find  a more  efficient  numerical  integration  procedure,  or  (2) 
modify  the  equations  in  such  a way  that  the  high-frequency  solutions 
are  eliminated.  The  first  alternative  leads  to  integration  methods  of 
implicit  type.  The  nonlinearity  of  the  resulting  implicit  system  and 
the  difficulty  of  incorporating  the  near  balance  of  the  equations  have 
not  been  successfully  dealt  with  so  far.  A strict  version  of  the  second 
alternative  is  actually  even  more  difficult  to  construct,  if  we  try  to 
eliminate  only  the  solution  of  sound-wave  type  and  obtain  a system  which  is 
still  hyperbolic.  The  main  difficulty  stems  from  the  fact  that  the  eigen- 
vectors corresponding  to  the  large  eigenvalues  are  different  for  the 
different  A.'s.  The  nonlinearity  of  the  system  and  the  effects  of 
variable  coefficients  are  further  complications. 

Various  approximate  versions  of  the  Eulerian  equations  have 
been  derived  by  intuition,  scale  analysis,  energy  conservation 
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considerations,  and  experience.  The  near  balance  of  the  equations  simplifies 
this  approach  considerably.  We  shall  study  two  such  approximate  sets  of 
equations:  the  hydrostatic  equations  and  the  incompressible  anelastic  system 
of  equations.  Unfortunately,  in  both  cases, the  hyperbolic  character  of 
the  system  is  lost. 

The  hydrostatic  system  is  derived  from  the  Eulerian  equations 
by  neglecting  the  vertical  acceleration  terms  in  the  third  equation  of 
motion.  One  so  obtains  the  "hydrostatic  equation" 

(3-1)  a I?  + g = 0 ' 

Here,  and  from  now  on,  we  use  the  notation  z for  as  a distinguished 

vertical  coordinate  and  w for  u^  as  the  vertical  velocity.  This 
approximation  is  extremely  accurate  for  the  large-scale  motion  of  the 
atmosphere.  The  wide-spread  use  of  the  hydrostatic  approximation 
actually  led  meteorologists  to  calling  the  resulting  system  "the 
primitive  equations  of  motion."  This  was  motivated  by  a comparison 
with  the  still  more  approximate  "quasi-geostrophic"  system,  but  the  term 
"primitive"  is  certainly  misleading. 

The  use  of  the  hydrostatic  approximation  has  several  important 
consequences.  First,  we  no  longer  have  a prognostic  equation  for  the 
vertical  velocity.  Second,  to  maintain  the  hydrostatic  equilibrium, 
the  time-changes  of  a and  p must  be  coordinated  in  such  a way  that 
a(dp/dz)  is  constant.  This  means  that  the  pressure  at  any  point  in 
the  model  atmosphere  can  be  determined  from  the  pressure  at  any  reference 
level  and  the  mass  of  the  separating  layer,  the  integral  / a dz. 
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With  this  approximation,  we  obtain  the  hydrostatic  system, 


system  Bl, 


(St  + *'v)%  + a VHP  f Eh  = 0 


a = o 


(^r  + jrV)a  - aV'u  = 0 
(J^-  + + pf7*u  = o 


where  u^Cu^Ug )* , u = (u^w^,  = (F^Fg)^  and  VR  = (d/dx^  d/dXg)1 

This  is  a much  more  complicated  system  than  it  may  seem  at  first  glance. 
The  equations  (3.2)  are  not  a hyperbolic  system.  To  show  this  we  form 
the  variational  equations 


+ + of?Hp'  + + «'VHp  + F^  = 0 


(3.3) 


-a'  f a' 
oz  a 


- o^7* ji'  + u' *va  - a'v*u  = 0 


+ u*v)p'  + prv-u'  + u' -Vp  + p'rv-u  = 0 


The  corresponding  constant  coefficient  problem  has  periodic  eigensolutions 
U1  = *1  exP^^vt  + + ^x2  + °^z ^ 3 > etc.,  for  large  and 

if  and  only  if 


(3.4) 


det  D ~ 0 


where 


-cd'-^lg+c  i 


prav  + ap 

„ xi  l 


-cd“2  (g  + Ciai^) 


pTa  + op 
*2  2 


dz  U1 


d 

2 2 
mz  ¥ c oP3 


praz  + «PZ 


cdA>  + — p ia>_ 
2 g x?  3 

p 

-a  (i+r)i“,v'u 
3 


- iv(g  + e ia>  ) 
6 ^ 


where  v = v + ^ ,u.u>.  + wo, . This  determinant  is  a cubic  polynomial  in 
J=1  J 0 3 

the  modified  eigenfrequency  V.  It  is  easy  to  see  that  the  roots  V , 

J 

j = 1,2,3,  of  (3*4)  have  the  asymptotic  behavior 


vx  = ©(1) 


v<:i)  • <X)  ■ • ■ - 


as 


Thus,  the  system  admits  solutions  with  arbitrarily 


large  signal  speeds  and  is  not  hyperbolic.  We  cannot  use  the  general 
methods  and  results  for  hyperbolic  systems  to  find  a well-posed  set  of 
boundary  conditions. 

Davies  [5  ] tried  to  avoid  this  unfortunate  effect  of  the  hydro- 
static approximation  using  a direct  energy  method  approach  similar  to 
the  one  used  in  Section  2 for  the  Eulerian  equations.  The  absence  of 
the  terms  coming  from  dw/dt  changes  (2.6)  to 
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(3.5)  ^ (a'n  RaH' 


■ ■ A 4 


- £ <s,*IV>  * (t  * 4 <EV*4  li^))q 


* »'*<£  (ba.1  - bc  - c*R>a'  - a'V'  - f'V 

where  g^  = (u^,u^,0,a' > V ) and  the  matrices  and  G are  obtained 

by  deleting  the  terms  arising  from  w(dw/dz)  in  the  original  third 
equation  from  the  matrices  and  C,  respectively.  Integrating 

(3*5)  over  the  region  fi  we  obtain 


(3-6)  ii-ll„.||2 


2 dt 


Hr 


< Kill'll2  * iia' Hr  11  - f 5^  (sjX^i  b JL 

+ (qj  *RA.j£'  ) dx^^dXgdz  . 


Davies  tried  to  find  a set  of  boundary  conditions  such  that  the  contri- 
bution from  the  boundary  integral  is  non-positive, proceeding  in  the  same 
way  as  we  did  with  the  Eulerian  system  A.  He  first  conjectured  that  the 
number  of  boundary  conditions  could  be  chosen  equal  to  the  number  which 
are  required  for  the  Eulerian  equations.  However,  this  conjecture  is 
false.  If  this  number  of  boundary  conditions  is  used  the  solution  cannot 
be  expected  to  satisfy  the  hydrostatic  relation  at  the  boundary.  The 


23 


r 


problem  is  overspecified  and  the  existence  of  continuous  solutions  is 
precluded.  There  is  a further  problem  with  this  approach.  The  energy- 
method  is  based  upon  norm  equivalences  and  the  norms  of  g'  and 
are  not  equivalent.  To  obtain  an  adequate  energy  estimate,  a bound  for 
the  term  g'  (RA^ ) - RC  - C R)q'  of  (3-5)  in  terms  of  ||g^|| 
instead  of  ||g'  ||  is  necessary,  but  this  is  not  possible.  Consequently, 
we  cannot  draw  any  conclusions  about  the  well-posedness  of  these  equations 
from  the  reduced  energy  equation  (3-6)  for  either  the  initial  boundary 
value  problem  or  the  Cauchy  problem. 

Since  the  energy  method  does  not  work,  we  now  turn  to  the  normal 
mode  analysis  technique. 

Normal  mode  analysis  of  the  general  equations  here  is  rather  complicated. 
If  the  motion  is  essentially  horizontal  we  may  instead  consider  linear- 
ization about  an  underlying  basic  state  a(z),  p(z)  which  satisfies 

op  + g = 0.  Due  to  this  simplification,  we  cannot  establish  sufficient 
z 

conditions  for  well-posedness  in  our  succeeding  analysis.  However,  we 
will  at  least  be  able  to  establish  some  necessary  conditions  and  we  can 
also  expect  the  variational  equations  to  reflect  the  main  properties 
of  the  system.  We  use  the  notation  a = a(z)  + a',  P + p(z)  + p'  and 
write  the  horizontal  velocities  as  u^  = v + where  v = (v^Vg)^ 

is  constant.  If  we  neglect  all  of  the  nonlinear  terms  in  primed  quantities, 
we  obtain  the  approximate  system 
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& »&  * ,Z,  Tj  ST  «h  * * £'  * 0 


s !£■-«*-<> 


it  °'  * vj  57  - S'7-*’  * V'  ’ 0 


Jf p’  + z vj  sir p* + prV’~  * ^ = 0 

J=1  J a 


We  can  then  transform  (3.7)  to  obtain  the  following  equations  in 
u^,  u^,  and  op', 


+ z vj  + vH(5p,)  + r = 0 


3=1  0 ~ 3 


(E  + Z v3  & l(«p'  } + VS,H"  0 


3=1  " j- 


where 


i<£®')  - -Is  - <OT’^  = - % - <&'» 

g \ pra  -g  £>z  / bz  a bz 


-1-2  ^ 

and  <J  = -g  a ^ In  9 is  the  "static  stability  of  the  basic  state. 

We  will  always  assume  a >0,  i.e.,  we  are  only  studying  perturbations 
of  a stable  stratification.  The  boundary  condition  at  z = 0 is 
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w'  = 0 which  implies  pra'  + op'  = 0.  As  an  upper  "boundary  condition 
we  have  p -*  0 as  z -»  ».  Our  condition  at  z = 0 can  also  be 
written  as  (op')  = gop'/a  so,  for  arbitrary  q' , 


f a_1q'  L(ap')  dz  = - / q*  C ' 

0 0 dz  \g2a 


dz 


a dq1  d 

g2a  dz  ^ 


(op* )dz 


(4 

ga 


q'^p* )z=0 


which  shows  th&t  the  operator  L is  self-adjoint  and  half-bounded. 
Therefore,  the  system  (3-8)  is  separable.  If  we  expand  the  variables 
u^  and  op'  in  the  eigenfunctions  F (z)  of  L,  we  then  obtain, 
for  each  V, 

+ £ TJ  SJ:1  »'k(v)  * V^tv)  + E(v)  * 0 

(3-9) 


V^v)  ” VH  • «H(V) 


+ V • u ' =0 


where  the  are  the  eigenvalues  of  L,  i.e., 


(3.10) 


LFv(z)  * Vv(z)  * 
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Since  o > 0,  this  eigenvalue  problem  is  of  Sturm-Liouville  type 


and  the  eigenvalues  are  distinct,  positive,  real  numbers.  For 

each  value  of  v the  system  (3.9)  is  hyperbolic  and  has  the  same 

form  as  the  shallow-water  equations  to  be  discussed  later.  The 

characteristic  velocities  are  v.,  v.  + c and  v,  - c where 

j j v j v 

cv  - \ • Under  standard  atmospheric  conditions  cQ  322  m/s., 

c^  53  3b  m/s.,  c^  33  17  m/s.  and  cv  -»  0 as  v ->  <*  , see  Wiin-Nielsen 
[22].  From  this  it  follows  that,  for  small  V,  two  of  the  character- 
istic velocities  will  be  positive  if  v.  > 0 and  one  will  be  positive 
if  v.  < 0.  When  v is  so  large  that  c < |v  |,  then  all  three 

J v j 

characteristic  velocities  are  positive  if  Vj  > 0 and  negative 
if  Vj  < 0. 


It  follows  that  the  appropriate  number  of  boundary  conditions 
for  our  simplified  problems  must  be  different  for  these  two  classes 
of  eigensolutions.  For  those  components  with  cv  < |v  |,  all  variables 
should  be  prescribed  if  vn  < 0 (inflow)  and  no  variables  should  be  prescribed 

if  vn  > 0 (outflow).  For  the  other  components,  the  appropriate  number  of 

boundary  conditions  is  two  if  vn  < 0 and  one  if  v^  > 0.  Possible 

forms  for  these  boundary  conditions  are  examined  later  for  the 

equivalent  shallow-water  equations.  The  simplest  choice  of  inflow 

conditions  is  probably  to  give  both  velocity  components , 

u^  , for  all  V and  c*P'(v)  for  those  values  of  v with 

c„  < |v  I.  At  outflow,  V > 0,  the  normal  velocity  component  can  be 
v n 7 n 

given  for  those  v with  cv  > |vj|.  In  the  special  case  of  a solid- 

wall  boundary  there  is  only  one  class  of  eigensolutions,  since  cv  > !v  I 

for  all  v.  The  condition  u (z)  =0  yields  u’  =0  for  all  V 

n n(v) 
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which,  as  shown  later,  yields  a well-posed  problem  for  each  v and, 
consequently,  for  the  entire  system.  We  reiterate  that  these  con- 
clusions are  only  valid  for  our  simplified  version  of  the  variational 
equations  and  that  we  have  no  proof  of  their  validity  for  the  complete 
system.  However,  the  conclusion  that  the  boundary  conditions  mist 
be  separated  in  terms  of  the  vertical  eigenfunctions  is  valid  for 
the  complete  system  since  that  system  admits  particular  solutions  of 
the  type  that  we  have  discussed.  Our  inability  to  provide  sufficient 
conditions  for  the  complete  system  is  essentially  due  to  the  fact 
that  we  cannot  show  that  all  solutions  of  the  complete  equations 
can  be  expressed  in  terms  of  the  eigenfunctions  of  (3.12).  It  follows 
from  this  discussion  that  local,  pointwise  boundary  conditions  cannot 
yield  a well-posed  problem  for  the  open  boundary  problem  for  the 
hydrostatic  equations;  well-posed  problems  can  only  be  obtained 
(l)  if  the  boundary  conditions  are  formulated  in  terms  of  local 
eigenfunction  expansions  or  (2)  nonlocal  boundary  operators  are  used. 
We  know  of  no  successful  formulation  of  the  second  type. 

A convenient  byproduct  of  the  hydrostatic  approximation  is 
the  possibility  of  using  variables  other  than  z as  the  vertical 
coordinate  (e.g.,  pressure,  p,  potential  temperature,  etc.).  If 
pressure  is  used  as  the  vertical  coordinate  the  equations  (3.2)  become 


(3-11) 


d_ 

dt 


Ct 


+ ^ = 0 

PT 


+ 

dP 


a = 


0 
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where 


— = A.  + u.  *V  + CD 
dt  dt  ^ p dp 


“ = ft  ' * = gZ  ’ 


and  V denotes  the  horizontal  gradient/divergence  operator  on 
constant  pressure  surfaces. 

The  p-system  has  the  advantage  that  the  region  of  integration 

has  a limited  vertical  extent,  0 < p < p , where  p is  the  surface 

s s 

pressure,  instead  of  0 < z < °°  for  the  original  system.  The 
meteorological  data  are  also  collected  and  analyzed  as 
functions  of  pressure  which  simplifies  the  construction  of  initial  data 
fields  on  constant  pressure  surfaces.  One  disadvantage  is  that  the  lower 
boundary  condition,  w = 0 at  z = 0,  becomes  dqp/dt  = 0 at  the  unknown 
surface  p = p (x,t)  where  <p  = 0.  The  usual  way  to  overcome  this 
difficulty  is  to  prescribe 


(3.1 2! 


if  “ < & + “ s’*  * 0 


at  a constant  pressure  surface  p = p^,  usually  chosen  as  lOOOmb, 

instead  of  at  p = p . The  upper  boundary  condition  simply  becomes 

s 


It  is  easy  to  show  that  the  transformation  to  the  p-system 
does  not  change  the  nonhyperbolic  character  of  the  equations,  and 
again  an  upper  bound  cannot  be  found  for  the  rate  of  growth  of 
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disturbances  on  the  solution  in  any  conventional  energy  norm.  As 
before,  we  have  to  limit  the  detailed  study  to  an  approximate  system 
of  variational  equations.  These  are  given  by 


<!r  - 4 * V f r 


= 0 


(3-13) 


(4  + E vi  dx“)a'  + (S  + - ) 

d,  J=1  j dx. 


dp  pr 


4 <p'  + a'  - 0 


¥ 


where  a = a(p).  We  are  only  considering  the  mean  translatory  part 
of  the  advection  terms.  The  system  (3-13)  will  retain  the  essential 
features  of  the  complete  variational  equations  if  the  solution  is  close 
to  steady-state  and  the  motion  is  quasi -horizontal. 

If  we  eliminate  to'  and  a' > the  system  (3. 13)  can  be 
written  in  terms  of  ul  and  T*  as  follows: 


(3.1k) 


z siHsi  + V + 1 - 


V5t  ' " * j dx  % p' 

J J 


(4  «■  £ v.  4-)L<p'  + V -u*  = 0 
ot  J=1  J d:tj  P "■« 
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I 


where 


JL  (ct_1  <*£'_) 

dp  ^ dp  ; 


+ jL) 

(5p  pt}  ' 


To  show  that  (3.14)  is  also  a separable  system  we  must  show  that  the 
upper  and  lower  boundary  conditions  are  consistent.  At  p = 0 we 
have  cu»  = o so  ^(dcp'/dp)  = 0.  At  p = p^  the  condition  (3. 12) 
may  be  transformed, using  the  equation  (3-13 ) for  a1 , to  obtain 


(3.15) 


so  that 


J J 


a-1  + 5"%'  =0  at  p = pQ 


if  this  condition  is  satisfied  initially.  Then 


f q'l/p'dp  = / 0 ^ dp  + (q,  "^q* cp* 


which  shows  that  L is  self-adjoint  and  half -bounded.  We  can  now 
expand  and  qp*  in  the  eigenfunctions  of  L.  The  simplified 

variational  equations  then  become  a family  of  hyperbolic  systems  , 


(•^r-  + Z v.  t-— )u,'  + v7  •p'  , + F!  \ 

^ j=l  J <*J  \v)  P (v)  (v) 


= 0 


(3.16) 


(■^“  + X)  V.  TT“~  )K  Cp|  \ + V U/  \ — 0 

dt  , x j dx . vMv)  p~(v) 


where  the  <v  are  the  eigenvalues  of  L,  i.e., 


(3-17) 


LGv(P)  = \Gv(p) 


with  corresponding  eigenfunctions  Gy(p). 

Again,  the  horizontal  boundary  conditions  must  be  formulated 

differently  as  the  normal  component  of  velocity  is  larger  or  smaller 
-1/2 

than  . All  of  the  conclusions  made  for  the  z -coordinate  system 

apply  in  this  case. 

We  now  consider  the  second  approximation  of  the  Eulerian 
equations  of  motion,  the  incompressible,  anelastic  system  (B2 ) . 
Approximations  of  this  type  have  been  used  in  many  areas  of  theoretical 
and  applied  fluid  mechanics.  In  their  most  simple  form,  the  equations 
for  an  incompressible  fluid  are 


(3.18) 


t aQ^p  = o 


V-u  = 0 


where  aQ  is  a constant.  This  system  is  usually  derived  from  basic 
physical  considerations, but  it  can  also  be  obtained  from  the  complete 
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Eulerian  equations  by  letting  a approach  a limiting  constant  value,  a^- 


■ 


The  more  general  anelastic  approximation,  often  used  in  studies  of  convective 
systems  in  the  atmosphere,  is  based  upon  the  following  assumptions: 

1)  the  potential  temperature  of  the  fluid  is  nearly  constant  (we 
denote  this  constant  value  by  5 in  our  following  discussion); 

2)  the  pressure  deviates  only  slightly  from  a hydrostatic  stratifi- 
cation; and  3)  the  typical  horizontal  and  vertical  length  scales 
are  similar.  If  the  characteristic  length  scale  is  much  smaller  than 
the  "scale  height,"  c^/g,  we  obtain  the  system 


& % + 8V  - 0 


Q 2 dw  p _ 

_ w + g _ . g _ = 0 


_d 

dt 


e 

$ 


(3.19) 


£0  = o 


V’u  = 0 


where  9 is  the  deviation  of  the  potential  temperature  from  a basic 
state  S,  7T  = Cp((p/pQ)  n - (p/pQ)  ' ),  p denotes  the  isentropic 
pressure  profile  corresponding  to  9,  and 

dt  = 5t  + %'VH  * W 3z  • 

We  note  that,  expect  for  the  equation  for  9 , this  system  is  essentially 
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the  same  as  I3.IB).  A third  system  of  the  same  basic  form  is  the 
Boussinesq  system  which  is  often  used  in  oceanography. 

The  linearized  variational  equations  corresponding  to  the  system 

(*.19)  are,  deleting  the  small  zero-order  terms, 


0-20) 


Sv.-ir' 

ri 


= 0 


A.  w.  + 9 

dt 


dir1 

dz 


- g 


0* 

— = 0 
s 


A + u' .ve=  0 

dt  ~ 

V*u  • = 0 . 

/V 


The  existence  of  nontrivial  periodic  eigensolutions  of  the  form 


exo[ i (vt  + ^1x1  + 

O)  x + 
2 2 

(jo^z ) ] is 

equivalent 

to  the 

determinant 

““  / 

iv 

0 

0 

9ia> 

N 

0 

\ 

/ 

0 

iv 

0 

9io>2 

0 

\ 

(5.21)  detl 

0 

0 

iv 

?iu> 

5 

-&/s 

r 

V 

ae_ 

dxL 

d*2 

dP 

dz 

0 

<v 

iv 

/ 

Vi^ 

ia>2 

icn3 

0 

o-y 

where 

V = V + 

2 

E 

3=1 

u.^.  + 
0 0 

3b 


Equation  (3-2l)  is  third  order  in  V and  of  fifth  order  in  the  to.. 

J 

If  the  lower  order  terms  are  neglected,  then  the  resulting  equation  has 

v = -Z2  ,u.a'.  - wu,  as  a triple  root.  Thus,  the  anelastic  eauations 

j=l  3 j 3 

are  not  a hyperbolic  system,  but  their  eigensolutions  have,  to  highest 

order,  time  dependent  behavior  which  is  like  that  of  hyperbolic  systems. 

We  cannot  obtain  sufficient  conditions  for  well-posedness  of 

the  initial  boundary  value  problem  for  time-singular  systems  like 

(j.20)  using  the  normal  node  analysis  technique.  The  theoretical 

justification  is  lacking  at  present.  However,  rigorous  results  on 

the  necessary  form  of  the  boundary  conditions  can  be  obtained  since 

it  is  clear  that  pathological  solutions  can  be  constructed  via  the 

normal  node  technique  following  Agmon's  construction  [12]. 

Analysis  of  the  eigensolutions  of  (3»20)  shows  that  four 

boundary  conditions  must  be  given  at  inflow  parts  of  the  boundary  and 

that  one  condition  must  be  given  at  outflow  parts  of  the  boundary. 

Furthermore,  pathological  behavior  like  that  exhibited  by  solutions 

of  the  approximate  system  B1  is  not  present.  It  can  be  shown  using 

the  energy  method  that  the  physical  boundary  condition  u^  = 0 

for  a rigid  wall  be.  ’ary  yields  a well-posed  problem  with  an  energy 

norm  of  type  (uj2  + U.J2  + w'2  + a2?5'2)'1'^2. 

For  the  constant  coefficient  problems  (3.20),  L estimates 
for  the  well-posedness  can  be  obtained  directly  using  Fourier-Laplace 
transform  techniques  for  the  initial  boundary  value  problem  on  a quarter-space; 
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t > O5  > 0;  -00  < x^,  z < 00 . It  follows  that  these  problems  are 
well-posed  if  y,  and  9 are  given  at  inflow  and  if  is  given 

at  outflow.  However,  as  mentioned  above,  the  reduction  of  the  general 
problem  to  quarter-space  problems  and  the  variable  coefficient 
problems  to  constant  coefficient  problems  via  freezing  arguments 
is  not  covered  by  existing  theory. 


1+ . The  Barotropic  or  "Shallow -Water"  Equations  (System  C). 

A third  approximation  to  the  Eulerian  equations , the  shallow- 
water  equations  which  are  our  system  C,  may  be  written 


(U.l) 


5T  % * V * £ ' 0 
i 9 + <??„■%  . 0 


where 


_d_ 

dt 


y 


and  Ujj  and  cp  are  functions  of  time,  t,  and  the  horizontal  space 
coordinates,  x^  and  x^.  F represents  any  zero  order  or  forcing 
terms  such  as,  e.g.,  the  coriolis  forces.  The  variable  <p,  the  geopotential, 
is  always  positive.  In  most  geophysical  problem,  the  flow  is  subsonic 
so  that  cp  > u^  + Ug,  but  the  opposite  relation  holds  for  both 
supersonic  flow  and  for  sane  of  the  subsystems  derived  from  the 
hydrostatic  equations  (Bl)  by  separation  of  variables. 
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At  least  three  different  names  exist  for  system  C.  In 
meteorological  applications  it  is  usually  called  "the  primitive 
barotropic  equations."  In  oceanography  the  most  common  name  is 
"the  shallow-water  equations."  This  system  accurately  describes  wave 
motion  on  the  surface  of  a homogeneous  fluid  when  the  horizontal 
wave  length  is  much  longer  than  both  the  vertical  scale  of  motion 
and  the  depth  of  the  fluid.  The  vector  form  of  system  C is 


(U.2) 


where  £ = (u.,^,^)1, 


1/2 

The  eigenvalues  of  A.  are  u.,  u.  + c and  u.  - c with  c = <p 

J J J J 

They  are  all  real  and  have  distinct  eigenvectors.  The  symmetric  and 
positive  definite  transformation  matrix 


*-l 


simultaneously  symmetrizes  A^  and  Ag.  Thus,  the  system  C 

has  retained  (or  regained)  the  basic  property  of  the  original  Eulerian 

equations  (system  A)  of  being  a quasi-linear  system  of  hyperbolic 


equations. 
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As  before,  we  are  only  discussing  problems  with  smooth 
solutions.  The  basic  properties  of  the  system  may  then  be  found  from 
the  corresponding  linearized  variational  equations 


IU.3) 


— n' 

dt  3, 


2 
- L 
J=1 


V*’ 


Ca*  + F'  = 0 


where  j is  a solution  of  ( b.2 ),  and  the  small  disturbance  a'  (x^,x?,t 
may  be  generated  by  the  inhomogeneous  term  F'  or  be  caused  by  an 
initial  disturbance  g/ (x1?x  ,0).  The  matrix  C has  the  form 


The  matrix  C,  the  coefficient  matrices  A^(a)  and  and  the 

transformation  matrix  R(a)  are  slowly  varying  in  space  and  time. 

It  is  easy  to  show  the  well-posedness  of  the  Cauchy  problem  in  the 
2 

L -norm  by  the  classical  energy  method.  This  was  done  by  Elvius 
and  Sundstrbm  [9]  and  Davies  [5]-  They  also  showed  the  well-posedness 
of  the  initial-boundary  value  problem  for  some  possible  sets  of 
boundary  conditions  using  the  energy  method.  However,  as  pointed  out 
by  de  Rivas  [ 7 ],  Davies  overspecified  the  boundary  conditions  in  his 
paper.  For  this  reason  our  discussion  follows  that  of  Elvius  and 


Sundstrbm. 
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The  growth  equation  for  the  energy  norm  of  g'  is 

<kA)  SF  £ *'*  Rl'  4x2 


Since 


2 

flt  + 5x7  (RAj } - RC  - C*R  la'  dXjL  d* 

'1  J-  J 


- ^ a.'  RF'  + F'  Rg/  dx.dx  - fa'*  M q ' 
« 1 2 dfi  n 


^9 


The  first  integral  is  bounded  by 


) + tp'2  )dx 


if  we  can  neglect  the  contribution  from  the  inhomogeneous  term  and 

if  the  boundary  conditions  ensure  that  the  integral  f-  2.,RAn~'dS  i 0 * 

* 

The  growth  rate  of  f g/  Rgf  dx1  dxp  is  then  bounded  by  the  maximum 

SI  1 

value  of  the  quantity 


-1-  - ► 

dx^  dx2  J 


2X1/2 


(5*9) 


du. 


du„ 


dx,  dx^ 


1 


che  difference  between  the  deformation  and  divergence  of  the  basic  flow. 
Results  of  this  type  may  also  be  found  for  finite  disturbances  (which 
was  not  possible  for  the  Eulerian  equations).  Without  any  linearization* 
the  complete  variational  equations  are 


2 

+ >:  Alg+g’) 

j=l  J 


+ Cg*  + F’ 


0 . 


Using  the  complete  transformation  matrix  R(<j  + we  get 


it  [ a’ V 


= f (cp  +■  cp*  ) (u|2  + u’2)  *■  <p'2  dx^dXg 
= - r Au^.  + u2)[(qp  + ep»  )(u^?  + u-2)  + <p'2] 


+ (cp  + <p' 


+ r 2cpf 

n 


■ f q'  RF*  + F1  Eg*  dx  dx  - f q'  SA  q’ds 
fl  * ~ ~ 14  12dfl~  n~ 


= [ (ui®,2)  + air  (^cp,2)  ^1^2 


2 

= <6  u'<p'  ds  , 


we  again  obtain  a growth  estimate 


14.6) 


|r  r a’*  RclL  ^ Or 


i -1- 

\\^53t^  bx2J  J " SxL  J 


< max 


* / a'  Ra' 

if  F’  = 0 and  the  boundary  integral  f £f  HA  - u^cp’  ds  > 0. 

3 n 

Note  that  we  now  have  a bound  for  the  growth  rate  which  is  valid  for 
disturbances  of  arbitrary  size  and  only  involves  the  deformation  and 
divergence  of  the  undisturbed  flow.  The  solutions  £ and  £ + <j’ 
may  thus  be  any  pair  of  solutions  to  the  shallow-water  equations. 

p p p 

For  finite  a’,  the  integral  f £*  Rq'dx  dx  = f (cp+cp*  )(u'  + u'  ) + <p*  dx^dx 

a 1 2 a 12  12 

2 

is  no  longer  the  square  of  an  L -equivalent  norm  for  £* , but  as  long 
as  <p  + cp'  (the  thickness  of  the  fluid  layer)  is  strictly  positive  we 
can  apply  the  corresponding  Liapunov  theorems.  Usually,  |<p' | is 
so  much  smaller  than  |<p|  that  the  difference  between 

( f £'  R(£  + £*  )£*  dx  dx  )^2  and  the  energy  norm  ( f £'  R(£)£'  dx-dx  )^2 

a a 

is  negligible. 

All  these  estimates  are  valid  if  and  only  if  the  contribution 

from  the  boundary  is  strictly  nonpositive.  For  the  linearized 

r * 

variational  equations  we  had  to  require  l £f  RA  £ ds  > 0,  i.e., 

dfi  n 


pop 

f v (cp(u*  + u!c)  + cp*  ) + 2(pu'cp'ds 
n n -L  ” 


bn 


= f {u  cpu'2  + ^ (c  + u ) (u'  + <p*/c)2  - % (c-u  ) (u1  - qp'/c)2)  ds  > 0 

1 n ~ -i_o  n n O n n ' — 


bn 


n n 
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For  the  complete  variational  equations,  the  integral  is  instead 

<5  (u  + u [ (9  + tp' ) (u'2  + u]2 ) + 9' 2 ] + (29  + cp'  )u'  9*  ds 
bil  n n n n 

= f { (u  + u'  ) (9  + 9'  )uj2 
dfi  n n x 

g-  (9  + 9' ) ( (9  + j 9'  )/c  + un  +■  u ) (u^  + 9*/c  )2 ) 

- i (9  + 9' ) ((9  + ^ 9*  )/c  - u - u*  ) (u'  - 9*/c  )2 ) ds 
2 2 n nn 

l/p 

where  now  c = (9  + 9* ) ' . We  again  study  the  two  different  types 
of  boundary  conditions  separately. 

1.  Solid-wall  boundaries.  Since  the  normal  velocity  vanishes  at 

the  wall  we  have  ] £*  RAn  £'ds  = 0,  and  the  matrix  has  only 

one  negative  eigenvalue  at  the  boundary.  The  well-posedness 

of  the  problem  then  follows  directly  (also  without  linearization) 
for  the  boundary  condition  ufi  = 0. 

2.  Open  boundaries.  From  the  number  of  negative  eigenvalues  of 

it  follows  that  at  an  inflow  boundary  (u^  < 0),  two  boundary  conditions 

should  be  given  in  the  subsonic  region  where  |unl  < c 

(|u  + u' | < (9  + 9* )/c),  and  three  conditions  in  the  supersonic 

n n — 2 

region.  At  an  outflow  boundary,  one  boundary  condition  is  required 
if  the  flow  is  subsonic,  and  no  condition  should  be  given  if  the 
flow  is  supersonic. 
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If  the  flow  is  supersonic,  the  boundary  conditions  should  be 
such  that,  for  the  variational  equations,  the  values  of  all  three  of 


the  variables  u^,  u^  and  cp'  are  prescribed  at  inflow  points. 

For  subsonic  flow,  the  inflow  boundary  conditions  should 

determine  the  value  of  two  quantities  of  type  u^  - a^(u^  + <P' /c) 

and  u'  - qp'/c  - a„(u'  + cp* /c ) since  u{  and  u’  -qp'/c  are  the 
n 2 n -L  n 

characteristic  combinations  corresponding  to  the  negative  eigenvalues 

of  Ar.  Not  all  such  combinations  give  well-posed  problems,  but  if 
2 2 

c + u + 2u  a,  - (c  - u )aZ  > 0,  the  limited  growth  rate  for  the 
n n 1 n 2 — 

energy  norm  of  q'  gives  a direct  proof  of  well-posedness . We  may, 

e.g.,  choose  ai  = ag  = °*  This  condition  can  be  achieved  by  pre- 

1 /? 

scribing  uj_  and  un  - 2qp  at  these  inflow  points,  see  Elvius 

and  Sundstrom  [9],  since  then  u^  = 0 and  u^  - 2 (cp  cp'  )^2  + 2cp^2 
= u^  - 2qp'/ (<p1//2  + (cp  + qp'  )^2)  = u^  - qp^/qp1/2  = 0.  Actually,  this 
boundary  condition  ensures  that  only  bounded  growth  can  result  for 
all  finite  disturbances,  as  long  as  min{cp1//2,  (qp  + <p* ) + u^ 
is  positive.  Inserting  the  complete  expressions  for  u^  and  u^ 
in  the  boundary  integral,  the  integrand  becomes 

(u  + u'  )[  (qp  + qp* ) (u'2  + u]2)  + qp'2]  + (2qp  + qp’  )u'qp' 

= <p* 2 { (un  + uM[l  + 4(<p  + qp'  ) ( (qp  + qp'  )^2  + cp1/2)  2] 

+ 2(2qp  > qp* ) ((qp  + qp'  W2  + qp1^2)"1) 

= <p*2  ( (<p  + qp*  )1/2  + qp1//2 ) "2  { (un  + uO(6cp  + 5cp'  + 2 (cp  + qp*  )1/2  cp1//2) 

+ 2 (2qp  + qp' ) ( (qp  + qp' )1/^2  + qp1//2^}  . 
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If  0 < qp  < qp  + cp' , the  value  of 


l/? 

2 (2cp  + cp*  ) ( (qp + 

6cp  + 5<P'  + 2 (qp  +-  qp*  )1'2 


1/2 


1/2 

is  always  larger  than  qp  while  if  0 < qp  cp'  < qp,  this  expression 
has  (cp  + cp'  )1/2  as  its  lower  hound.  For  u^  +•  u^  + Min  (qp^2,  (qp  +■  cp*  )^2 ) > 0, 
the  inflow  part  of  the  boundary  then  provides  a nonnegative  contri- 
bution to  the  boundary  integral. 

Another  possible  inflow  boundary  condition  is  to  prescribe 

both  u and  ui,  so  that  a,  = 0,  a„=  -1.  This  condition  was 
n -1-  12 

apparently  first  suggested  by  Rousseau  [l6J.  The  value  of  la^l  is 
then  so  large  that  we  cannot  use  the  classical  energy  method,  but 
as  first  shown  by  Elvius  and  Kreiss  (private  communication),  well- 
posedness  can  be  proved  using  Kreiss*  normal  mode  analysis  technique. 

An  alternative  boundary  condition,  discussed  in  the  paper  by 
Elvius  and  Sundstrom  [9]  is  to  give  u_j_  and  cp  at  inflow  points.  This 
corresponds  to  the  choice  a = 0,  a =1.  For  problems  in  only 

one  space  dimension  this  is  a well-posed  condition.  The  value  of  ap 
is,  however,  so  large  that  we  cannot  use  the  energy  method  as  above. 

For  the  complete  two-dimensional  problem,  Elvius  and  Sundstrom  did 
not  analyze  the  well-posedness  properties,  but  numerical  experiments 
indicated  that  it  might  actually  be  an  ill-posed  set  of  boundary 
conditions.  This  conjecture  has  been  confirmed  by  a complete  analysis 
by  Elvius  and  Kreiss  (private  communication). 
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At  outflow  parts  of  the  boundary  (ufi  > 0),  we  should,  in 

the  subsonic  case,  give  one  boundary  condition.  For  the  variational 

equations  this  condition  should  prescribe  the  value  of  a combination 

of  the  form  u^  - cp’/c  - b^  - b2(u^  + <p'/c).  As  before,  well- 

posedness  follows  directly  from  the  positivity  of  the  integrand 
* 

$i'  RAnSL'  i-n  the  growth  equation  if 

(c  + un^un  + f'/c)^  + 2unui.2  ~ (c  ~ un^biuj_  + b2^Un  + <P,/C^2  > 0 

for  all  u*  and  u^  + <p* /c.  The  simplest  choices  of  the  parameters 
b^,  b^  satisfying  this  condition  are:  l)  b^  = b£  = 0,  i.e.,  giving 
un  ” <P'/C>  may  be  achieved  by  specifying  un  - 2qp  ; ; 2)  b^  = 0, 

bg  = -1,  i.e.,  giving  u^;  and  3)  bL  = 0,  bg  = 1,  i.e.,  specifying  <p 
also  yields  a well-posed  problem. 

All  of  these  conditions  guarantee  bounded  growth  for  finite 
disturbances.  In  each  case,  the  integrand  <j'  + 

strictly  positive  as  long  as  u^  + u^  > 0. 

For  the  open  boundary  problems  we  have  several  possible 
sets  of  boundary  conditions  which  all  satisfy  the  necessary  and  sufficient 
well-posedness  conditions.  If  the  problem  is  part  of  a telescoping 
technique  or  nested  integration,  or  if  any  arbitrary  type  of  boundary 
data  can  be  obtained  from  measurements , the  choice  between  these 
different  possibilities  may  be  difficult. 
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The  experiments  reported  by  Elvius  and  Sundstrom  [9]  do 

not  show  any  large  differences  between  the  results  from  a numerical 

model,  using  either  boundary  conditions  giving  the  value  of  the  combination 

u'  - ro’/c  and  u!  at  inflow  and  u'  - cp'/c  at  outflow  or  u'  and  u' 
n -L  n ^ n -L 

at  inflow  and  u'  at  outflow.  Their  results  examined  a long-wave 
n 

solution  with  small  long -wave  or  short-wave  disturbances. 

Further  experiments  by  Elvius  (private  communication)  indicate 
that  for  solutions  with  a less  pronounced  long-wave  character,  the 
first  alternative  is  less  susceptible  to  boundary  disturbances  which 

may  arise  when  the  normal  velocity  is  small  and  changes  sign. 

In  one  of  the  first  papers  on  limited-area  integration  of  the 

shallow-water  equations,  Charney  [2]  suggested  a quite  different 
set  of  boundary  conditions.  Since  it  is  both  inefficient  and 
difficult  to  implement,  this  set  is  now  primarily  of  historical 
interest.  It  is  still  worth  analyzing,  since  it  illustrates  the 
hazards  of  intuitive  deductions.  At  that  time,  most  experimental 
and  routine  work  on  numerical  weather  prediction  was  done  with 
"balanced"  forecast  models.  The  simplest  version,  the  barotropic 
vorticity  equation,  can  be  considered  as  a further  simplification  of 
the  shallow -water  equations.  In  the  derivation  of  this  approximation, 
system  C is  first  transformed  into  a set  of  three  prognostic 
equations  for  the  divergence  D = V*u,  the  vorticity 

£ = du2/dx^  - du^/dx^,  and  cp,  respectively  by  differentiating  the 
equations  of  motion  with  respect  to  x and  y and  combining  the 
results.  This  differentiated  system  is  then  simplified  by  using  a 
steady-state  approximation  in  the  first  and  third  of  these  equations 
and  by  keeping  only  the  nondivergent  advection  terms  in  the  vorticity 
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equation.  By  this  approximation,  the  whole  system  is  condensed  into 
one  prognostic  equation  in  one  dependent  variable , the  stream  function 
For  this  vorticity  equation,  Charney,  Fj^rtoft,  and  von  Neumann  [3] 
concluded  that  two  boundary  conditions  should  be  given  at  inflow  points 
of  an  open  boundary  and  one  condition  at  outflow  points.  They  suggested 
that  the  stream  function  (and  thus  the  normal  component  of  the  velocity) 
should  be  specified  at  all  boundary  points,  and  in  addition,  the 
vorticity  at  inflow  parts  of  the  boundary.  It  is  easy  to  show  that 
their  conclusion  on  the  number  of  boundary  conditions  was  correct  and 
that  the  suggested  conditions  make  the  problem  well-posed,  cf. 

Sundstrom  [19]. 

Charney' s proposed  boundary  conditions  for  the  shallow-water 
equations  were  apparently  based  on  the  idea  that  since  the  number  of 
boundary  conditions  is  the  same  for  the  vorticity  equation  and  the 
shallow  water  equations,  the  type  of  conditions  should  be  similar. 

He  therefore  suggested  that  un  should  be  prescribed  at  all  boundary 
points , and  as  the  second  quantity  to  be  given  at  inflow  parts  of  the 
boundary  he  chose  the  "potential  vorticity" 

du2  duL 

dx  " dx_  + f 

p = — xr~ • 

<P 

If  the  differentiated  version  of  the  shallow-water  equations  is 
formulated  in  terms  of  P,  D,  and  cp,  one  of  the  three  equations  is 
(d/dt )P  = 0.  This  differentiated  system  requires  three  boundary 
conditions  at  inflow  points,  one  more  than  the  number  required  for  the 
single  vorticity  equation,  and  one  of  these  conditions  may  be  the 
prescription  of  P.  However,  this  is  not  a valid  argument  for  the 
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usefulness  of  the  potential  vorticity  as  a boundary  condition  for  the 
undifferentiated  system-  We  cannot  just  pick  two  of  the  three 
necessary  inflow  conditions  for  the  differentiated  system  and  expect 
them  to  form  an  appropriate  set  of  boundary  conditions  for  the 
shallow-water  equations. 

The  danger  of  using  the  boundary  conditions  that  Charney 
suggested  can  be  shown  directly.  For  the  variational  equations,  they 
give  u^  = 0 and  (d/dn)uj_  - Pep'  =0  if  un  < 0.  The  last  condition 
cannot  be  used  directly,  but  it  can  be  combined  with  the  prognostic 
equation  for  u^  to  yield  an  equation  for  the  inflow  boundary  values  of 
(d/dt)u^  which  involves  only  boundary  quantities.  The  tangential 
velocity  is  then  determined  by  integrating  this  equation  from  t = 0. 

This  is  a very  complicated  way  of  computing  the  inflow  values  of 
Uj_.  Additionally,  this  approach  has  the  liability  that  a small  error 
committed  initially,  or  at  any  later  time  t^,  will  influence  the 
boundary  values  at  all  later  times.  These  errors  will  spread  into  the 
region  of  integration  contaminating  the  solution. 

5 . Effects  of  Viscous  Terms 

As  described  in  the  introduction,  small  viscous  terms  are 
often  added  to  the  systems  of  equations  we  are  considering.  They 
are  often  introduced  to  provide  a dissipative  filter  for  a numerical 
approximation.  In  other  cases  there  is  a physical  motivation  for 
using  viscous  terms  to  represent  diffusive  transport  (eddy  flux) 
of  momentum  and  heat.  Both  the  viscosity  and  heat  conduction  coefficients 
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are  usually  very  small,  but  since  these  terms  change  the  character 
of  the  differential  equations,  we  have  to  reinvestigate  the  boundary 
conditions . 

For  system  A (the  Eulerian  equations)  the  general  viscous 
form  of  (2.l)  is  the  compressible  Navier-Stok.es  equations,  here 
written  in  the  special  form 

— u +-  Q^p  + F = atu'^u  + (A  + uMv-u)] 

^ a - aV-u  = KnP^VCap1^) 

it  P + = 0 

where  u and  >,  are  the  Lame^^  constants,  and  where  the  heat  exchange 
term  represents  eddy  flux  of  potential  temperature.  Using  the 
potential  temperature  0 = (p^~  n /?,)  op  ' 1 as  dependent  variable 
instead  of  a,  the  second  equation  may  then  be  simplified  to 


V*9 


Here,  we  shall  only  study  the  simple  viscous  form  obtained  when 
= -A  = Va*^.  In  vector  notation,  we  have 


with  a = (u1,u2,\i5,0,p)t, 
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This  is  an  incompletely  parabolic  system,  see  Strikwerda  [18].  As 
before,  the  matrices  A^,  A^,  A^,  and  B can  be  simultaneously 
symmetrized  by  multiplication  from  the  left  by 
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where  a is  the  same  positive  parameter  which  occurred  in  the 
transformation  matrix  R in  section  2.  If  the  solution  (and  then 
also  R and  the  coefficient  matrices)  are  slowly  varying  in  space 
and  time,  we  can  use  the  general  results  of  Strikwerda  [18]  from 
which  it  follows  that,  for  well-posedness , it  is  necessary  and 
sufficient  that  the  family  of  linearized  variational  equations 


+ £ A (a) 

j=l  J 


0Xj 


fora  a well-posed  problem. 

We  first  study  the  rigid  wall  problem.  A normal-mode 
analysis  shows  that  four  boundary  conditions  are  required.  One  of 
these  conditions,  u^  = 0,  follows  immediately  from  the  solid-wall 
condition  ur  = 0.  The  three  remaining  conditions  may  be  chosen  as 

aiuii + (1-ai^  It  uii  = °>  a2ui2  * (1-a2)u  It  ul2  = °> and 


a^0'+  (l-a^)K^  ^ = 0 for  some  nonnegative  a^,  a^,  and  a^. 
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If  = 1,  these  conditions  correspond  to  a "nonslip", 

perfectly  conducting  wall,  while  if  a.^  = = a^  = 0,  they  represent 

a "perfect  slip",  thermally  insulated  wall.  The  well-posedness  of 
the  solid-wall  problem  with  these  boundary  conditions  can  be  demonstrated 
by  the  energy  method. 

For  small  values  of  u and  K^,  we  can  expect  boundary  layers 

of  thickness  & (li  ^ ) and  ) at  the  solid  wall.  If  some 

n 

a^  / 0,  the  value  of  the  corresponding  variable  can  be  expected  to 
change  by  a finite  amount  within  the  boundary  layer.  In  such  cases, 
we  must  use  a numerical  approximation  that  resolves  the  boundary 
layer  satisfactorily.  If  all  a^  = 0,  the  variables  will  only  change 
by  an  amount  proportional  to  the  boundary  layer  thickness,  and  the 
difference  between  the  solutions  to  the  viscous  and  inviscid 
equations  is  then  always  small  for  small  u and  K^. 

In  problems  with  open  boundaries , there  are  no  obvious 
physical  boundary  conditions.  As  before,  we  have  some  freedom  in 
choosing  the  most  suitable  mathematical  boundary  conditions.  For  the 
viscous  equations,  it  is  not  sufficient  to  prescribe  a set  of  con- 
ditions that  makes  the  problem  well-posed.  Since  we  know  that  there 
should  be  no  boundary  layers  at  the  inflow  and  outflow  boundaries , 
we  have  to  choose  the  mathematical  boundary  conditions  accordingly. 

This  problem  has  been  studied  by  Gustafsson  and  Sundstrom  [11]. 

They  showed  that  for  u,  -»  0,  the  conditions  must  yield  a well- 
posed  set  of  boundary  conditions  for  the  inviscid  equations,  i.e., 
singular  boundary  layers  should  not  occur. 
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appropriate  for  the  inviscid  equations.  Using  the  energy  method,  we  can  show 
that  the  following  set  of  boundary  conditions  gives  a well-posed  problem 
with  no  essential  boundary  layer: 

i / du' 

u^  - p'(a/pr)1/2  = o,  u|L  = u|9  = <?’  = 0,  and  u = = 0. 


12 


A similar  set  of  conditions  consistent  with  the  other  type  of  inviscid 

conditions  (u1  = u!,  = u’  = q'  = 0)  is  less  trivial  to  find.  This 
n JL.1 

is  because  none  of  these  inviscid  conditions  include  the  pressure 
disturbance  p’ . 

At  an  outflow  boundary,  the  viscous  equations  require  four 
boundary  conditions  and  the  inviscid  problem  only  one.  If  we  prefer 
to  give  u^  = 0 as  the  inviscid  condition,  we  can  now  prescribe 

"i  - °>  * S uIl  ■ **  ST  \ S 9'  “ 0 

The  inviscid  condition  u^  - p' (a/py)  7 =0  can  also  be  extended 

to  the  viscous  case  by  giving 


u;  - P^(a /pr)l/2  + n(a/pr)l/2  = °>  n ulr  n |r  = h 


0*  = 0. 

on 


In  both  these  cases,  we  can  use  the  energy  method  to  prove  well- 
posedness,  and  the  absence  of  boundary  layers  as  n,  ^ 0 

follows  directly. 
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All  these  boundary  conditions  may  be  expressed  in  terms  of 
a'  and  p'  instead  of  0'  by  using  the  identity 


e»  = e 


(°L  f 

Va  p r / 


which  is  valid  for  small  perturbations. 

For  the  hydrostatic  approximation  (Bl)  to  these  Euler ian 
equations,  the  analysis  is  at  least  not  simplified  by  including 
viscous  terms.  Even  if  we  could  find  a set  of  boundary  conditions 
that  make  the  linearized  variational  equations  with  constant  coefficients 
well-posed,  we  would  not  know  if  this  is,  in  some  sense,  true  for 
the  complete  equations. 

For  the  anelastic  equations  (B2),  viscous  terms  may  easily 

be  included.  To  find  the  proper  number  of  boundary  conditions,  we 

may  as  before  study  only  the  linearized  variational  equations  with 

constant  coefficients.  The  result  of  the  analysis  is  that  at  all 

boundaries,  four  boundary  conditions  should  be  prescribed.  At  an 

inflow  boundary  this  is  the  same  number  as  required  by  the  inviscid 

equations.  We  may  use  the  same  conditions  as  before,  u*  = u.'  = u!  = 0'  = 0 

n 11  ±2 

for  the  linearized,  constant  coefficient  equations. 

At  an  outflow  boundary,  the  inviscid  equations  required  one  boundary 

condition,  u’  = 0.  This  may  now  be  supplemented  by 
n 


P dn  Uil  4 dn  U12  ^ dn  ° * 
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For  the  solid-wall  case,  we  may  either  choose  the  inflow  or  outflow 
type  of  conditions.  The  choice  depends  on  whether  we  have  nonslip 
or  perfect  slip  and  a perfectly  conducting  or  a thermally  insulated 
wall.  The  well-posedness  of  chess  conditions  is  easily  demonstrated. 

For  the  shallow-water  equations  (c),  we  have  the  same  type 
of  behavior  as  with  the  Eulerian  equations  (A).  At  a £olid  wall, 
the  obvious  condition  u^  = 0 must  now  be  supplemented  by  one  more 
boundary  condition,  aiui  + (l-a^)u  U1  = ^1:ie  e^ec^s  caused 

by  choosing  a^  = 1 and  a^  = 0 are  similar  to  the  results  for 
system  (A). 

At  an  inflow  boundary,  the  boundary  condition  u^  - <p'/c  = 0, 

uj  = 0 is  easily  modified  by  adding  the  required  third  condition 

u u*  = 0.  If  the  inviscid  conditions  are  u’  = 0,  u'  = 0,  we 

* on  n n JL 

cannot  find  a suitable  viscous  form  by  our  energy  method  analysis. 

At  outflow  boundaries,  we  need  two  conditions  which  can  be 

chosen  to  be  u'  = 0 and  u u'  = 0 or 
n ^ on 

, du'  ^ 

u*  - SL  + M.  a = o and  u u'  = 0 , 
n c c on  ^ on 

depending  on  what  type  of  condition  is  preferred  for  the  inviscid 
equations.  Well-posedness  is  easily  proven  by  the  energy  method, 
and  no  artificial  boundary  layers  are  generated. 


I 
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Implications  for  numerical  methods 


The  analysis  of  various  initial-boundary  value  problems  contained 
in  this  paper  was  motivated  by  difficulties  arising  in  computational 
models  of  the  problems  we  have  discussed.  It  is  appropriate  to  comment 
on  the  implications  these  results  have  on  these  numerical  models. 

The  stability  analysis  and  related  error  estimates  for  the 
approximate  methods  are  the  discrete  analogues  of  our  well-posedness 
analysis  and  estimates  of  the  form  (l.j)*  For  a given  well-posed  problem 
of  the  types  discussed  here,  we  can  always  find  stable  difference 
approximations  and  numerical  boundary  conditions.  Examples  are  given 
in  Gustafsson,  Kreiss  and  SundstrSm  [10]  and  Elvius  and  Sundstrom  [9]. 
Conversely,  an  approximation  cannot  have  a norm  which  behaves  reasonably 
if  it  accurately  approximates  an  ill-posed  problem. 

When,  for  a given  problem,  the  number  of  boundary  conditions  is 
overspecified,  the  difference  approximation  may  well  be  stable.  However, 
the  effective  boundary  conditions  which  influence  the  solution  are, 
in  general,  difficult  to  determine,  especially  for  problems  in  several 
space  dimensions.  They  may  well  be  a complicated  function  of  the 
conditions  given  and  bear  little  resemblance  to  them. 

An  additional  complication  induced  by  overspecification  is  that  the 
underlying  solution  being  approximated  is  not  generally  continuous.  The 
phenomena  associated  with  approximations  to  discontinuous  solutions  have 
been  studied  by  several  authors,  a good  discussion  and  summary  of  these 
results  can  be  found  in  section  10  of  Thomee  [21].  These  results  may 
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be  summarized  as  follows.  If  a non -dissipative  approximation  is  used, 
then  high  frequency  waves  emanate  from  the  region  of  the  discontinuity 
and  travel  across  the  domain  without  losing  appreciable  amplitude. 

They  will  usually  travel  with  the  highest  fundamental  wave  speeds  and 
rapidly  cover  the  domain  with  error.  If  scalar  equations  are  being 
approximated,  then  this  region  of  error  can  be  restricted  to  the  vicinity 
of  the  discontinuity  by  using  dissipative  approximations.  However,  these 
results  do  not  apply  to  systems  of  equations  as  we  have  here.  The  errors  can 
propagate  away  from  the  discontinuity  through  other  components  of  the 
solution.  Boundary  value  overspecification  may  be  regarded  as  a 
stationary  source  of  such  discontinuities. 

In  order  to  avoid  the  problems  associated  with  the  proper  selection 
of  boundary  conditions,  the  order  and  type  of  the  differential  equations 
is  often  raised  to  obtain  a problem  that  is  easier  to  analyze  and 
approximate.  The  equations  are  usually  modified  by  adding  dissipative 
terms  so  that  the  number  of  boundary  conditions  is  appropriate.  Unfortunately, 
this  idea  seldom  works.  If  a spurious  boundary  layer  of  appreciable  size 
results,  the  effects  are  not  unlike  those  discussed  above  for  discon- 
tinuities and,  unless  the  dissipative  terms  are  very  large,  the  error 
introduced  at  the  boundary  will  again  propagate  into  the  interior. 

If  the  boundary  conditions  are  underspecified  there  are  no 
a priori  estimates  for  the  differential  equations.  In  order  for  an 
approximation  to  be  computable  there  must  be  a sufficient  number  of 
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boundary  conditions  specified  for  the  approximation.  This  cannot  be 
fewer  than  the  number  required  for  the  differential  equation. 

Additional  conditions  are  usually  constructed  by  means  of  extrapo- 
lations. For  an  underspecified  problem  the  extrapolation  of  quantities 
that  should  be  prescribed  results  in  an  unstable  or  inconsistent 
approximation  of  the  wrong  differential  equation. 
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